Filamentation processes and dynamical excitation of light condensates in optical 

media with competing nonlinearities 
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We analyze both theoretically and by means of numerical simulations the phenomena of filamen- 
tation and dynamical formation of self-guided nonlinear waves in media featuring competing cubic 
and quint ic nonlinearities. We provide a theoretical description of recent experiments in terms of 
a linear stability analysis supported with simulations, showing the possibility of the observation of 
modulational instability suppression of intense light pulses travelling across such nonlinear media. 
We also show a novel mechanism of indirect excitation of light condensates by means of coalescence 
processes of nonlinear coherent structures produced by managed filamentation of high power laser 
beams. 
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I. INTRODUCTION 

Studies on nonlinear beam filamentation of ultrashort 
pulses are being actively developed nowadays |l| . Far be- 
yond a certain critical power threshold, extremely short 
pulses may become unstable and break up into several 
uncorrelated light structures by modulational instabil- 
ity (MI) processes. These structures, called filaments in 
this context, are spatiotemporal soliton-like light distri- 
butions which appear after MI and hold for large propa- 
gation distances without changing their shape and size[2j. 
In fact, optical filaments have been observed in several 
scenarios like air[l], water [3], carbon disulfide [J], etc., 
and some mechanisms to control their spatial location 
have been recently proposed [5]. 

On the other hand, the theoretical description of 
these localized modes is in general based on nonlinear 
Schrodinger equations (NLSE)[6|-[8j. In the modeling of 
these kind of systems, nonlinear effects like two and three 
photon absorption, second order group-velocity disper- 
sion, time-delayed nonlinear response, multiphoton ion- 
ization and plasma defocusing[9] are usually taken into 
account. Nevertheless, in order to favor the insight of 
dynamics in these systems, simpler models are employed 
when certain physical conditions are satisfied. To this 
end it is common to consider NLSE involving just lo- 
cal intensity-dependent nonlinear susceptibilities which 
are real[||. As an example, the so-called cubic-quintic 
(CQ) nonlinear ity[10| has been used as a solvent approx- 
imation in modelling subpicosecond pulse propagation in 
nonlinear media [4]. In the latter reference, the dynamics 
of filament formation and evolution in carbon-disulfide 
cells is analyzed both experimentally and by numerical 
means. At this respect, extensive numerical work based 
on a conservative CQ model is successfully employed in 
the qualitative description of experimental achievements. 

In addition "liquid light condensates" , i.e, robust soli- 
tonic distributions in CQ optical media with a "flat-top" 
transverse spatial envelope and intriguing surface tension 
properties were studied in (ll) . Recently, it was shown 



that their dynamics follows the same equations governing 
the evolution of usual liquid droplets jl2|. These results 
were obtained by both analytical and numerical meth- 
ods, but they have not yet been confirmed by real exper- 
iments. The main difficulty, of course, is the search for 
a cubic-quintic nonlinear medium, although some candi- 
dates were inspected in |13j- In fact, both nonlinear third 
and fifth order terms in the polarization of the material 
are usually complex, so it is often necessary to consider 
some other nonlinear processes like ones related above. 
Nevertheless, in the light of Ref.[4] we are allowed to 
think of materials like CS2 as real CQ nonlinear me- 
dia in specific power regimes. Recently, some common 
gases like air, N2 and O2 have also been identified as 
CQ optical media by means of the measurement of 77,4 
coefficients [li). although the propagation dynamics of 
probe pulses throughout these media has not yet been 
investigated. 

The purpose of the present paper is threefold. Firstly, 
we will give a simple theoretical description of the results 
obtained in Q by means of a linear stability analysis of 
plane waves. Secondly, we will study both analytically 
and numerically the possibility of real observation of the 
modulational instability suppression process in CS2 • Fi- 
nally, we will provide a detailed description of a novel 
mechanism of dynamical excitation of flat-top solitons in 
"pure" CQ systems via filamentation management and 
coalescence and we will analyze the possibility of achiev- 
ing light condensates in carbon disulfide. 



II. PHYSICAL MODEL 

We will consider the propagation of a high-intensity 
laser pulse through a CS2 bulk, within the framework of 
the nonlinear Schrodinger equation (NLSE). In order to 
simplify the theoretical model, we assume a scalar slowly 
varying spatial envelope and we neglect group-velocity 
dispersion effects. These approximations have been ex- 
perimentally justified for media like carbon disulfide [4] 
and air [9] within proper parameter spaces. We will also 
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include a cubic-quintic nonlinearity, resulting in the fol- 
lowing NLSE 

2ik n ^ + ViV> + 2k 2 n (n 2 \^\ 2 - n 4 |^| 4 )^ = (1) 



III. LINEAR STABILITY ANALYSIS OF PLANE 
WAVES 

It is well known that Eqj2] admits plane wave (PW) 
solutions of the type $ = $oe* 7r7 , where &o is a constant 
which depends on 7 eigenvalues as follows: 



where the quantities above are defined as follows: 
ip is the complex slowly varying electric field envelope 
of the corresponding electromagnetic wave propagating 
through the nonlinear medium, ko = 27r/A is the vac- 
uum wavenumber corresponding to the laser source wave- 
length A. no is the linear refractive index of the medium 
and z is the propagation distance. Nonlinear coefficients 
722,724 characterize the strength of the real third and 
fifth order nonlinear optical susceptibilities respectively. 

= d 2 /dx 2j rd 2 jdy 2 is the Cartesian transverse Lapla- 
cian operator. Notice that EqJI] preserves the norm of 
the field i.e. we have also neglected the effect of non- 
linear two-photon absorption, which can be justified for 
moderate intensity pulses travelling throughout a carbon 
disulfide bulk over a short distance [4[ . 

For technical reasons related with numerical treatment 
of partial differential equations, it is often useful to ex- 
press Eq.(pQ) in a reduced form by introducing a suitable 
set of non-dimensional variables, namely 



$0 



1 

2 2 



(3) 



Let us proceed with a linear stability analysis of these 
homogeneous solutions, following the Bespalov-Talanov 
procedure . The first step of this method is to add a 
small amplitude perturbation of the form £ = £r + z£j 
to the PW solutions, with £ RiI = f0^ e [iK X T-HM and 
\£,R,i\ <C $0, where Ki is the transverse wavevector and 
h is the propagation constant of the corresponding per- 
turbational mode. Now, for h 2 < 0, depending on the 
sign of we have two possible outcomes. If 5s (h) < 0, 
the perturbations grow exponentially in 77 yielding to field 
destabilization. However, if $s(h) > the amplitude of 
the perturbational mode decays quickly, so that the un- 
derlying PW remains stable. On the contrary, whenever 
h 2 > 0, the amplitude of the perturbations remains lim- 
ited. 
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where the adimensional variables considered from now 
on are: 77 = zkon^/n^ is the adimensional propagation 
distance, A = d 2 /dx 2 + d 2 /d( 2 is the dimensionless 
transverse Laplaci an opera tor, being the spatial coordi- 
nates O = 0koU2 y^noAu, with O = \-> C? = x^y. The 
normalized beam irradiance |<£| 2 = (714/722) 1^1 2 is mea- 
sured in units of the quotient 724/722 which univocally 
characterizes the nonlinear response of the medium. 

Keeping in mind the equivalence between EqsJT] and 
[2j typical parameters corresponding to beam propaga- 
tion in real dielectric materials like CS2 are (for fem- 
tosecond optical pulses with A = 800nm|15|) no = 1.6, 
n 2 = 3 • 10- 15 cm 2 /W. To our best knowledge, the CS 2 
quintic coefficient 724 has not been measured yet, since 
there are no values of this quantity cited in the literature, 
although an approximation 724 = 2 • 10 ~ 27 cm 4 /W 2 based 
on the agreement between simulations and experiments 
was successfully introduced in [4(. Very remarkably, the 
first experimental determination of the nonlinear coeffi- 
cient 724 of several gaseous media has been recently car- 
ried out in [lij . 

In the following section, we will introduce an analyt- 
ical linear stability analysis of plane waves in the CQ 
model and its results will be compared with those of the 
experiments in Q. 
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FIG. 1: Plot of Ymax vs. |<3>o | 2 . The highest growth rate cor- 
responds to <I>o = 0.25, whereas for $0 > 0.5 all values Ymax 
vanish or become negative delimiting the domain of modu- 
lational instability suppression. The labeled points on the 
graph refer to the pseudocolor intensity plots on the right side 
of the picture. These plots display a small region [width=500, 
height=300] of perturbed PW of intensities $0 = 0-045 (a), 
$0 — 0.25 (b), $0 — 0.64 (c) at the propagation distance 
77 = 125. The filaments shown in (b) arise before those in (a), 
so that both coalescence and domain formation effects can be 
observed. 

By substituting the expression for the perturbed PW 
$ = ($ + i)e %171 in Eq|2]and separating both the real 
and imaginary parts of the outgoing relation, we obtain, 
after linearizing on £r,£/, i.e. retaining only the first 
order terms in perturbation theory, the following couple 
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of Schrodinger-type equations 



1.25 



^-A^-(2*g-4#g)6 



drj 



■ = 



(4a) 



(4b) 



where £r 5 j depend upon K± and ft, in the way noted 

above. From Eqj4b] we obtain £r = and taking 

into account Eqj4a|we deduce the algebraic relation 



T = K 1 _[-Ki + 2$l 



■4*3] 



(5) 



where V is the growth rate of the perturbations, satis- 
fying r = ih. Notice that all the information about the 
PW perturbational spectrum is contained in the ampli- 
tude <3>o- We can estimate the value of the transverse 
wavevector which gives the maximum growth rate by 
looking for the extremals of Eq. [5j i.e 



dT 

dK ± 



(K™ ax f = §1 



(6) 



As a consequence, the PW of amplitude $o = \ fea- 
tures the highest perturbation growth rate (see FigJT]). 
This means that it will suffer destabilization in a shorter 
characteristic scale than any other homogeneous solu- 
tions of Eq. [2j if it is perturbed by a harmonic mode 
with K± = K™ ax . These assertions are illustrated by 
FigsJT]|2] In FigJT] we have plotted the quantity T max 
as a function of the PW intensity. For |$o| 2 > 0, we 
get T max < 0, corresponding with the modulationally 
stable PW branch [17j. The three pseudocolor intensity 
plots in the right column snapshots of FigJT] display dif- 
ferent examples of the state of perturbed light homoge- 
neous distributions after a propagation distance n = 125. 
All PW considered throughout this manuscript are per- 
turbed with a 5% spatial random noise and propagated 
with a standard beam propagation method. 

The top picture FigJT]- (a) corresponds to a simulation 
0.045 (unstable PW). It is observed that this 
noisy distribution leads to the formation of several self- 
trapped optical channels which survive a long distance. 
We have chosen this intensity because it corresponds to 
one of those used in the experiments We estimate 
the FWHM of the arising filaments to be around 10/im, 
in good agreement with the experimental measurements. 
The middle picture FigJT]- (b) shows the last stage of a 
simulation with $q = 0.25, corresponding to the PW 
with the highest T max . In fact, as this PW posseses the 
shortest longitudinal scale of destabilization, it can be 
seen in FigJT]- (b) that after filament formation some do- 
mains are formed due to coalescence processes between 
self-guiding structures. The final example in the bottom 
picture FigJT]- (c) shows a stable PW of $q = 0.64. This is 
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FIG. 2: Plots of the peak intensity evolution over the line 
[C = 0,-120 < x < 120] of perturbed plane waves. Blue-solid 
line [<J>o — 0.64] is almost constant, featuring small bounded 
oscillations due to spatial noise initially added. Orange- 
dashed line [3>o — 0.25] shows a rapid growing in intensity 
with 77, which is a trace of a short longitudinal instability 
scale. Black-dotted line [$o = 0.045] shows the destabiliza- 
tion scale for experimental data[4]. It is important to notice 
the difference in between both instability scales of orange- 
black curves. For the simulation with experimental data, the 
starting point of the filamentation process is around 77 = 75 
(2.1mm in dimensional units) which is compatible with that 
of experiments. In all cases, the maximum propagation dis- 
tance is n = 300. 



a very interesting illustration of the phenomenon called 
modulational instability suppression (MIS). For $q > 0.5, 
the linear perturbational modes decrease exponentially, 
thus leading to a suppression of the modulational insta- 
bility which could be observed in real experiments. To 
our knowledge, this would constitute the first example 
of beam stabilization in a nonlinear optical medium ob- 
tained just by changing the input intensity. 

Moreover, when the considered perturbations are in- 
tense enough, the quantity T max acquires a greater rele- 
vancy since its inverse value yields the characteristic lon- 
gitudinal scale of instability development A||[16j. In a 
first approximation, we will consider the filamentation 
regime to be described by the coefficient A|| as well as 
by the quantity A^ = 7r/K™ ax , which univocally char- 
acterizes the transverse spatial scale of the emerging 
filaments 0, [l6j|. We are assuming throughout this pa- 
per that optical filaments are axially symmetric so that 
A± is valid for both transverse coordinates (x, C). I n 
general, both quantities (A||, A_l) are useful in estimat- 
ing the order of magnitude of different spatial instability 
scales, but their agreement with both simulations and ex- 
periments can be rather qualitative due to the existence 
of different filamentation regimes, radial perturbations, 
nonlinear responses, inhomogeneities in the medium, coa- 
lescence processes between filaments, etc. In any case, as 
multiple filamentation appears as a consequence of mod- 
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ulational instability of the corresponding field, the linear 
stability analysis should be at least qualitatively correct 
at the leading order, where other nonlinear processes are 
not considered to be involved. 

In Figj2j we have selected a transverse cut of the three 
paradigmatic PW considered above (see FigJT]) and we 
have determined their peak intensity along the propaga- 
tion. As it can be appreciated, the perturbed PW with 
<£q = 0.64 (blue-solid line) features a quasi-constant in- 
tensity during the whole simulation, illustrating the MIS 
process. However, for <3>o = 0.25 (orange-dashed) and 
$o = 0.045 (black-dotted), a dramatic growing of the 
peak intensity with 77 is observed, indicating the destabi- 
lization of the initial conditions. Furthermore, there is a 
substantial difference in the longitudinal instability scale 
of both situations. It can be appreciated that $0 = 0.25 
(PW featuring the highest T max ) displays shorter insta- 
bility scale than <3>o = 0-045. Hereafter, we will com- 
pare this numerical result with the theory as well as the 
predictions of the stability analysis for the experimen- 
tal parameters considered by Centurion et al. [J], that 
observed the filamentation process undergone by fem- 
tosecond pulses propagating in CS^- Using pulses of 
peak intensity I p = 0.68 GW/mm 2 , they measured a 
mean filament sizes of 12/im, in good agreement with 
our estimation. Moreover, the transverse scale of the fil- 
aments predicted by our linear stability analysis for their 
experimental configuration is A± = tt/(^q — 2<£q)2, or 
Aj_ « 16.45/ira in dimensional units, which is also con- 
sistent with the sizes measured in [4]. This can be in- 
terpreted as a further theoretical proof of the validity of 
the CQ model in describing the dynamics of femtosecond 
optical pulses propagating in CS2- 

On the other hand, the longitudinal scale of instability 
development A|| can be estimated as A|| « = 24.4 

(0.7mm in dimensional units). This scale has not been 
directly characterized by Centurion et al. in Ref. [4j. 
Nevertheless, they found completely developed filaments 
on a scale of 5mm, which is compatible with our value of 
A||, that constitutes a rough estimation of the propaga- 
tion distance for the filaments when they start to grow 
up. Note that for <I>q = 0.25 we get A|| =8, thus justify- 
ing the shorter characteristic scale of instability develop- 
ment observed in Fig 121 

While our analytical and numerical results for A^ are 
both in good agreement with the experiment of Ref. [4j , 
the estimations for A|| turn out to be less accurate. This 
is due to the fact that in the simulations we have per- 
turbed the homogeneous distributions with random noise 
instead of using harmonic perturbations with fixed K±. 
In fact, we have checked that using harmonic perturba- 
tions leads to smaller values for A|| by a factor of ~ 3, 
getting closer to the analytical calculation. 

These results show that, in spite of its simplicity, the 
linear stability analysis not only provides a good insight 
into the physics, but it can also be used to obtain quan- 
titative useful predictions for the experiments of filamen- 
tation in CS2- 
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FIG. 3: Intensity pseudocolor plots of the SGF evolution in 
the presence of an unstable homogeneous background. The 
intensity snapshots correspond to a spatial domain [—75 < 
X, C < 75] at several propagation distances, namely: (a) 77 = 0; 
(b) j] = 80; (c) j] = 150; (d) 77 = 200; (e) 77 = 1000. Last pic- 
ture (f) displays the phase map of the field depicted in (e). 
From this phase structure, it can be stated that the light con- 
densate appearing at the last stage of the propagation is a co- 
herent light structure since the color within the region where 
it is located (see blue circle in picture (f)) is homogeneous. 



IV. LIQUID LIGHT SOLITONS ARISING BY 
COALESCENCE 

In this section we will show through numerical simu- 
lations a novel mechanism for the dynamical excitation 
of so-called "liquid light condensates" [llj in media with 
CQ nonlinear ity, by means of spatial control of the fila- 
mentation regime and coalescence processes between self- 
guided optical channels. Furthermore, we will discuss 
about the possibility of observation of these structures in 
the experiments with CS^- 

As starting point, we have considered an initial condi- 
tion for our simulations consisting on a broad-extended 
laser pulse with random noise fluctuations added. Group- 
velocity dispersion effects shall not be considered here 
since it is justified to assume that the temporal profile 
of the pulse does not change during short propagation 
distances through CS2 [II • We will model this light dis- 
tribution with a wide fiat-top profile represented by the 
following function; 



n = A 



bg 



0.25 ( 1 + tanh(p + c^) 



1 



tcmh(p-uj b g) ) 
(7) 

where A^g is the peak amplitude of the light pulse, 

9 — Vx 2 + C 2 is the radial coordinate and Ub g is the 
mean radius. We have performed all simulations fixing 
A^g = 0.2 and uo^g — 1000. As we will restrict ourselves 
to a small region located at the center of the pulse, we 
can think of the domain delimited by Q as part of a plane 
wave, since boundary effects (e.g. near-field diffraction) 



will not affect our results within the longitudinal scale 
imposed by the common experimental setups. We will 
explore the case A^g = 0.2, corresponding to a modu- 
lationally unstable plane wave [17] as we have discussed 
above, although the same qualitative behavior is obtained 
whenever T max > (see FigJTJ. 

By adding an inhomogeneity to the light distribution 
of Eq. (|7|), we will be able to manage the fllamentation 
regime. In fact, we have induced a small scale Gaus- 
sian fluctuation (SGF) of amplitude Aq — 0.5071 and 
width uoq = 30, located at the center of the bidimen- 
sional spatial profile Q in order to select a small region 
where fllamentation will first appear. Notice the scale 
difference between the two counterparts. Then by over- 
lapping both structures ft and SGF, we are generating a 
spatially inhomogeneous profile &i n with a peak intensity 
|^in| 2 = 0.5, which corresponds to that of the plane- wave 
defining the boundary between the stable and unstable 
PW solutions of Eqffl 0. 

The evolution dynamics of such light distribution is il- 
lustrated in Figj3j In this figure, the snapshots show our 
small region of interest [—75 < xX < 75]. The initial 
condition is plotted in picture (a). Once the pulse has 
entered the medium, the SGF starts to self-focus. It is 
remarkable how such a process is fastly counteracted by 
the self-defocusing nonlinearity. As a result, a bright hot 
spot is observed at the distribution centroid and diffrac- 
tive rings representing Shockwaves due to "background- 
fluctuation" interplay appear (see (b) snapshot). 

The spectral properties of these waves (usually referred 
as conical emissions) have been studied in [18[. They give 
support to the formation of new filaments from pulse MI. 
Hereafter, we will discuss how to control the shockwave 
propagation in order to stimulate the appearance of sev- 
eral filaments surrounding the hot spot (see (b) snap- 
shot). The key point is to generate a SGF in a precise 
way by using phase masks [5] and techniques like pulse 
shaping [l9j]. To excite wide liquid solitons, we propose 
to surround the main custom-generated light structure, 
which is roughly similar to a top-flat beam, with a crit- 
ical number of filaments located close of it. In fact, the 
solitons are statistical attractors due to their high in- 
ternal energy. The satellite filaments are then used to 
provide an energy source for the central seed so that its 
power could be increased by coalescence processes, if the 
phases of the soliton and the filament match. As shown in 
snapshot (c) in Figj3j there is indeed a dynamical energy 
exchange between filaments and hot spot during propaga- 
tion. As a result, all solitonic structures are both created 
and annihilated several times (picture (c)), until a wider 
coherent structure with almost constant peak intensity 
arises (picture (d)). As a consequence, the emerging liq- 
uid light soliton [llj is strongly perturbated by the energy 
excess released in the collective coalescence. Neverthe- 
less, from snapshot (e) we see that after a large propaga- 
tion distance the flat-topped soliton features an improved 
radial symmetry, the energy excess being radiated in the 
form of linear modes. 
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FIG. 4: Transverse intensity profile of the filaments distribu- 
tion (black solid line) at n — 5000, with a superposed trans- 
verse profile of an eigenstate (red dashed line) of Eq. [2] with 
propagation constant 7 = 0.173. It is seen that the agree- 
ment between both contours is reasonable. FWHM (delim- 
ited with blue solid lines) of the different-sized filaments rep- 
resented in the figure have been estimated to be L c « 28.5 
(flat-topped soliton) and Lf « 10 (Kerr- type soliton). The 
spatial coordinate % G [—150, 150]. 



In picture (f) we also show the phase map correspond- 
ing to the intensity plot displayed in (e). We see that 
within the region where the condensate is placed (high- 
lighted with a blue circle) the phase structure is homo- 
geneous. 

Finally, we have checked the transverse size of the dis- 
tribution after a huge propagation distance, n = 5000. 
This allows the light condensate to relax towards a sta- 
tionary state, thus minimizing random profile variations 
due to surface modes. In Fig|4j we have compared this 
extracted profile with that of a stationary solution of 
Eq|2j The agreement between the transverse profiles cor- 
responding to this simulation and the stationary mode is 
remarkable. This completes the demonstration that this 
nonlinear structure observed in our simulations is a true 
flat-topped CQ-soliton. 

To determine whether this kind of top-flatted trans- 
verse modes have been observed in the experiments with 
CS^^y], we have estimated their transverse size, and 
found that it would lie around the value L c & 28.5 
(30.2/im in dimensional units). This is almost 3 times 
greater than the size of the optical filaments measured in 
Ref. [4]. We thus conclude that Centurion et al. have 
not reached the threshold to produce top-flatted, liquid 
light solitons. 

On the other hand, on the right side of the light con- 
densate displayed in Fig 0] we also see a completely devel- 
oped self-guided light channel of size Lf « 10 (10.6/ira 
in dimensional units), which resembles closely to those 
observed by Centurion et al. This is a further demon- 
stration that our simulations reasonably reproduce the 
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existing observations, besides providing a guide for the 
new managed production of the higher power top-flatted 
solitons. 



V. CONCLUSIONS 

In this paper, we have provided theoretical support to 
recent experiments on filamentation of ultrashort pulses 
propagating through carbon disulfide by means of ap- 
proximate analytical methods and numerical simulations 
of propagation. The agreement between our calculations 
and the experimental outcomes is remarkable, taking into 
account that we have developed the theoretical descrip- 
tion assuming a quintic coefficient which has not been 
measured in the laboratory. We have shown that un- 
der certain conditions, a suppression of the modulational 
instability could be observed in such media over an inten- 
sity threshold. We have also described a novel procedure 
for indirect excitation of flat-topped solitons in pure CQ 
materials, and we have argued that these nonlinear struc- 
tures have not been observed experimentally yet. Very 



remarkably, recent experiments with air and other gases 
have revealed their CQ nature by measuring their high- 
order nonlinear coefficients [l^. Thus our results open 
the door to the quest for liquid light condensates in real 
experiments. 
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